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A dynamical many-body theory is presented which systematically extends beyond mean-field and perturba- 
tive quantum-field theoretical procedures. It allows us to study the dynamics of strongly interacting quantum- 
degenerate atomic gases. The non-perturbative approximation scheme is based on a systematic expansion of 
the two-particle irreducible effective action in powers of the inverse number of field components. This yields 
dynamic equations which contain direct scattering, memory and "off-shell" effects that are not captured by the 
Gross-Pitaevskii equation. This is relevant to account for the dynamics of, e.g., strongly interacting quantum 
gases atoms near a scattering resonance, or of one-dimensional Bose gases in the Tonks-Girardeau regime. We 
apply the theory to a homogeneous ultracold Bose gas in one spatial dimension. Considering the time evolution 
of an initial state far from equilibrium we show that it quickly evolves to a non-equilibrium quasistationary state 
and discuss the possibility to attribute an effective temperature to it. The approach to thermal equilibrium is 
found to be extremely slow. 
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I. INTRODUCTION 

Since the pioneering achievement of Bose-Einstein conden- 
sation in dilute alkali-metal gases much effort has been made 
to extend this success to a wider class of atomic and molec- 
ular species as well as to a wealth of different trapping ge- 
ometries. In particular, theproduction of quantum degenerate 
ensembles of molecules [1, 2] promises a wide spectrum of 
important applications ranging from ultra-precise molecular 
spectroscopy and cold collision studies yO to "superchemi- 
cal" reactions fT] and the investigation of the crossover from 
a Bose-Einstein condensate (BEC) of molecules to Bardeen- 
Cooper-Schrieffer (BCS) correlated pairs in Fermi gases jST. 
Furthermore, (quasi) one- and two-dimensional traps ji^Q] 
as well as optical lattices J3, 01 allow us to realize strongly 
correlated many-body states of atoms reminiscent of similar 
phenomena in condensed matter systems. 

Zero-energy scattering resonances,^ particularly the so 

called magnetic Feshbach resonances llOt [ill nA uM so 
far have played a leading role in the creation of strong inter- 
actions in degenerate atomic quantum gases. Near a Feshbach 
resonance, the scattering of, e.g., a pair of Bose-condensed 
atoms, whose relative energy is very close to zero, can be 
described by a strongly enhanced s-wave scattering length 
a. Present-day experimental techniques allow for resonance- 
enhanced scattering lengths larger than the mean interatomic 
distance n^^/'^ in the gas. As a consequence, the diluteness 
parameter = na!^ is larger than one. The Bose-Einstein con- 
densate is no longer in the coUisionless regime, it represents a 
strongly interacting system. 

Feshbach resonances also play an important role in the 
physics of ultracold, degenerate Fermi gases, where they al- 
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low US to induce a transition from a phase of Bose-Einstein- 
condensed diatomic molecules to a superfluid phase of BCS- 
type where the Fermion pairs are correlated over large dis- 
tances 1 15, il6l IITT . In the transition region, the interaction 
is strong and the correlation length easily exceeds the mean 
interatomic distance. 

A strongly interacting system can also be obtained with- 
out the requirement of a scattering resonance: In a one- 
dimensional trap, the gas enters the so-called Tonks-Girardeau 
regime, if the dimensionless interaction parameter 7 = 
gmm/ {h^n) is much larger than one [18, 19, 201. Here, giD 
is the coupling parameter of the one-dimensional gas, e.g., 
(71D = 2fi^a/{m.l'\) for a cylindrical trap with transverse 
harmonic oscillator length l^. In the Tonks-Girardeau limit 
7 ^ 00 the atoms can no longer pass each other and behave 
in many respects like a one-dimensional ideal Fermi gas |7]. 

In an optical lattice, strong effective interactions can be in- 
duced by suppressing the hopping between adjacent lattice 
sites and thus increasing the weight of the interaction rela- 
tive to the kinetic energy [8. .21.] . This leads, in the limit of 
near-zero hopping or strong interactions, to a Mott insulating 
state 122]. 

Conventionally, the dynamics of many-body systems of 
suffiently weakly interacting particles is described using per- 
turbative approximation schemes for the exact quantum field 
theoretical many-body equations of motion. Such schemes are 
in principle based on expansions in terms of powers of some 
dimensionless parameter, like y^, which measures the binary 
interaction strength. It is clear that for strongly interacting 
systems for which the relevant parameter is no longer small, a 
perturbative approach must eventually fail. 

A particular challenge represent systems far from equilib- 
rium. For instance, the sudden change of an external magnetic 
field for atomic quantum gases near a Feshbach resonance 
can lead to dramatic non-equilibrium phenomena i23l I24I1 . 
If the magnetic field is modified on time scales smaller than 
a typical collisional duration, the system is driven far from 
thermal equilibrium. Standard approaches based on small 
deviations from equilibrium, such as linear response theory. 
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or kinetic descriptions requiring a sufficient homogeneity in 
time, are not applicable. Another example is the formation 
of a Bose-Einstein condensate through controlled evapora- 
tive cooling of an ultracold atomic gas, an every-day exper- 
imental procedure, which has been studied in great detail 
(cf. e.g. Refs. a complete theoreti- 

cal description of which is, however, still lacking fsP]. 

In this article a dynamical many body theory is presented 
which systematically extends beyond mean-field and pertur- 
bative quantum-field theoretical approximation schemes. The 
non-perturbative approach is based on a systematic expansion 
of the two-particle irreducible (2PI) effective action in p ow- 
ers of the inverse number of field components I32L 133113411 . 
The 2PI 1/M expansion to next-to-leading order yields dy- 
namic equations which contain direct scattering, memory and 
"off-shell" effects. It allows to describe far-from-equilibrium 
dynamics as well as the late-time approach to quantum ther- 
mal equilibrium. 

Recently, these methods have allowed important progress 
in describing the dynamics of strongly interacting relativistic 
systems far from thermal equilibrium for bosonic 1 33, 35, 36, 
l37>.38.1 as well as fermionic degrees of freedom |39, 40]. Our 
aim is to employ the 2PI effective action for ultracold quan- 
tum gases and to numerically solve the 2PI expansion to 
next-to-leading order This is exemplified for a homogeneous 
ultracold Bose gas in one spatial dimension. We compute the 
time evolution of an initial state which is far from thermal 
equilibrium. After a characteristic short-time scale it is found 
to be driven to a quasistationary state. However, the system 
is still far from equilibrium and the thermal equilibration time 
can exceed the early-time scale by orders of magnitude. In 
particular, a unique temperature can not be attributed to the 
quasistationary state. This is important in view of experi- 
ments with one-dimensional traps, where the longest acces- 
sible times may be too short to see complete thermalization. 

We emphasize that mean-field approximations fail to de- 
scribe the dynamics even qualitatively. Similarly, standard ki- 
netic descriptions based on two-to-two collisions give a triv- 
ial (constant) dynamics because of phase space restrictions 
in one spatial dimension. To our knowledge this is the first 
time that the one-dimensional dynamics has been described 
from a full systematic 2PI expansion in a non-relativistic 
quantum field theory. The 2PI effective action approach in 
this context has been discussed previously using an additional 
weak-coupling expansion and solved for a Bose gas in a lat- 
tice in one spatial dimension |41]. See also Ref. |42] for a 
discussion of the 2PI (or ^-derivable) approach to the theory 
of weakly interacting Bose gases. The connection to kinetic 
theory (cf., e.g., Refs. 143*, 44, "4?, "46", "4?, 48, 49, 50, "sT] 
and references therein) has been discussed in Refs. 152, 53.1. 
cf. also Refs. 01 El El E3 El in in the context of ele- 
mentary particle physics and cosmology. The 2PI 1 /M ex- 
pansion has also been successfully applied to compute critical 
exponents in thermal equilibrium near the second-order phase 
transition of a model in the same universality class |60|. The 
good convergence properties of the expansion for small val- 
ues of M > 1 have also been observed in the context of non- 
equilibrium classical statistical field theories 16 LI . 



Ultracold atomic Bose gases in one- and two-dimensional 
trapping configurations, due to their specific quantum statisti- 
cal properties, have been studied in great detail, using both 
perturbative and non-perturbative approaches. In lower di- 
mensional systems, the regime of applicability of perturba- 
tion theory (cf., e.g., ll62l l63ll ) is reduced as compared to 
three-dimensional gases since fluctuations, in particular of the 
phase, play a stronger role such that non-perturbative descrip- 
tions are required. Stationary systems have been studied us- 
ing, e.g., renormalization group ll64ll65ll or Monte-Carlo tech- 
niques |66, 67]. 

Our article is organized as follows: In the remainder of this 
chapter we introduce the basic principles and ideas underly- 
ing the non-perturbative approximation scheme. In Chapterllll 
the 2PI effective action is defined, and the (exact) dynamic 
equations are deduced from it. We then show, in Chapter 
Unl that the well-known dynamical Hartree-Fock-Bogoliubov 
equations of motion result from a single diagram in a loop 
expansion of the 2PI effective action. In Chapter CVl we in- 
troduce in detail the non-perturbative approximation scheme 
on the basis of a 1/M expansion of the 2PI effective action, 
where M is the number of field components in a subspace 
where the action is invariant under 0{N) rotations. For the 
simplest case of a single complex Bose-field, we have M =2. 
Before we draw our conclusions in Chapter IVTl we apply, in 
Chapter |V] the many-body dynamic equations, to next-to- 
leading order in the 1/M expansion, to a uniform ultracold 
Bose gas in one spatial dimension. 

A. Non-perturbative approximations out of equilibrium 

The systems of identical bosons to be considered in this 
article are described by the Hamiltonian 

j d^x*^(x)iJiB(x)\['(x) 

+ ly d\d?y¥{^)¥{y)V{\^~y\)^{y)i>{^), 

(1) 

iJlB(x) = -— ^+ytrap(x), (2) 

The second term on the right hand side of Eq. Q represents 
the bare interaction term. The bare coupling constants equal 
the binary interaction potential T^(|x — y|) at the possible in- 
terparticle distances |x— y|. This potential is obtained, e.g., by 
describing the quantum mechanics of two-atom interactions in 
the Born-Oppenheimer approximation, where |x — y | denotes 
the internuclear distance. 

The dynamics of the non-relativistic many-body system can 
be described, e.g., in the Schrodinger picture, where the many- 
body state at time t is given by some density matrix poit)- 
The Schrodinger equation, with the Hamiltonian Q, then 
completely determines the non-equilibrium dynamics. All in- 
formation about the quantum theory can be encoded in the 
infinite series of correlation functions or ri-point functions 

{i!\xi) ■ ■ ■■^{Xn))t ^l:r[pD{t)¥ {Xl) ■ ■ - -^{Xn)], (3) 
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with n > 1. For the Hamihonian Q, the equation of motion 
of any particular rt-point function involves other correlation 
functions up to the order of n + 2. Hence, the result is an infi- 
nite system of coupled dynamic equations for the correlation 
functions (cf., e.g., Ref. 168. 691). Since this cannot be solved 
exactly one has to find suitable approximation schemes. 

Here it is important to note that for out-of-equilibrium cal- 
culations there are additional complications which do not ap- 
pear in vacuum or thermal equilibrium. The first new aspect 
concerns secularity: Even for weak couplings the strict per- 
turbative time evolution suffers from the presence of spuri- 
ous, so-called secular, terms which grow with time and in- 
validate the expansion. Moreover, the very same problem 
appears as well for non-perturbative approximation schemes 
such as standard 1/J\f expansions based on the one-particle 
irreducible (IPX) effective action |88|. Similar problems can 
also appear by simply truncating the infinite system of cou- 
pled dynamic equations for the correlation functions at a given 
level of n-point functions. 

The problem of secularity has been discussed in Refs. 
fyol |7l'|. Typically, for a given approximation, there can 
be various ways to resolve the secularity problem by resum- 
mation. There is a requirement, however, which poses very 
strong restrictions on the possible approximations: Universal- 
ity, i.e. the insensitivity of the late-time behavior to the details 
of the initial conditions. If thermal equilibrium is approached 
then the late-time result is universal in the sense that it be- 
comes uniquely determined by the conserved energy density 
and particle number. To implement the necessary nonlinear 
dynamics which recover detailed balance at late times is de- 
manding. Both requirements of a non-secular and universal 
behavior can indeed be fulfilled using systematic expansions 
of the 2PI effective action, which provide a practical means to 
describe far-from-equilibrium dynamics as well as thermal- 
ization from first principles rTOl lTlll . 

As will be described below, to lowest order, if all quantum- 
statistical fluctuations are neglected, the 2PI effective action 
leads to the Gross-Pitaevskii equation for the mean field fs^. 
For a sufficiently large mean field this classical field-theory 
approximation can be used to approximately describe the dy- 
namics of weakly interacting Bose-Einstein condensates. Of 
course, the classical approximation cannot be used to describe 
the approach to quantum thermal equilibrium characterized by 
a Bose-Einstein distribution. 

The 2PI effective action, with fluctuations taken into 
account to two-loop order, describes the (time-dependent) 
Hartree-Fock-Bogoliubov(HFB) approximation, in which ex- 
change between the condensate and the non-condensed frac- 
tion of the gas is accounted for. It conserves total particle 
number and energy but neglects multiple scattering. As a 
consequence, the approximation fails to describe thermaliza- 
tion. Typically, the maximum time for which it is reliable, de- 
creases with increasing interaction strength, e.g., in the case 
of ultracold alkali atoms, scattering length a. In particular, the 
Hartree-Fock approximation is known to suffer from the pres- 
ence of an infinite series of additional conserved quantities, 
which are not present in the fully interacting theory. These 
spurious constants of motion are associated to an infinite life- 



time of quasi-particle momentum modes, which prevent re- 
laxation to a thermal distribution |33|. 

The extensively employed HFB equations of motion may 
also be obtained by rewriting the hierarchy of dynamic equa- 
tions for the correlation functions (|3} in terms of their con- 
nected counterparts, or cumulants, and neglecting all cumu- 
lants of three and more field operators |'68', "t^I. One finds 
that the HFB dynamic equation for the leading order cu- 
mulant, the mean field, involves a resummation of infinitely 
many graphs and therefore arbitrary high powers of the bare 
coupling V. For Bose-Einstein condensates, many approxi- 
mation schemes beyond HFB have been introduced, e.g., in 
Refs. |43, 68], and take into account correlation functions of 
third and higher orders. However, the requirements for a re- 
liable late-time behaviour are difficult to implement. The ap- 
proximation schemes, first, have to yield controlled approx- 
imations and, second, must not violate crucial conservation 
laws like energy conservation, or, for a non-relativistic gas, 
the conservation of total particle number 

We emphasize that all these requirements can be met by 
systematic expansions of the 2PI effective action. These in- 
clude 2PI loop-expansions, coupling-expansions, as well as 
1/A/^-expansions. In the following we concentrate on the ex- 
pansion in inverse powers of J\f, where Af is the number of 
field components of a scalar theory invariant under the sym- 
metry transformations of the orthogonal group 0{Af). Com- 
pared to a coupling expansion, this procedure has the advan- 
tage that it can be employed in the absence of a weak coupling. 
In particular, it can be used to describe the non-analytic dy- 
namics near the second-order phase transition |60], where the 
condensate vanishes — a situation which cannot be quantita- 
tively described in terms of a coupling expansion. The scalar 
theory defined by Eq. Q is J7(l)-invariant, which is equiva- 
lent to an 0(2)-symmetry. In contrast to the standard ex- 
pansion of the IPI effective action, the apparently rapid con- 
vergence of the 2PI expansion even for small values of JV is 
crucial for our approach |60, 61]. In the following we present 
the 2PI effective action approach to the dynamics of a strongly 
interacting ultracold Bose-gas. 



II. 2PI EFFECTIVE-ACTION APPROACH 

A. Non-equilibrium quantum field theory 

Before introducing the 2PI effective action we would like 
to recall, for the purpose of making our article sufficiently 
self-contained, some basics about non-equilibrium quantum 
field theory. For more details, cf., e.g., Ref. ]71]. All infor- 
mation about a non-equilibrium quantum many-body system 
is contained in the generating functional for non-equilibrium 
correlation functions I9(J1: 



Z[J,K;pD]^Tr\pD{to)'rcexph f ^x)J{x) 

J X 



(4) 
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Any correlation function of a quantum many -body system can 
be derived from this by functional differentiation with respect 
to J{x) and subsequently setting J = K = 0. For example, 
the two-point function follows as 



5^Z[J, K; pd] 



i5J{x)iSJ{y) 



(5) 



.J=K=0 



Here, poito) is the normalized density matrix describing the 
many -body system at the initial time to, which, in general, 
does not have the equilibrium form p^'^'' ~ exp{— /3i/}. $ 
is the operator of a TV^-component scalar quantum field, with 
its space-time arguments x = {t, x) ~ (xq, x). In the follow- 
ing, where obvious, we suppress indices i enumerating the 
two components $i, corresponding to the real and imaginary 
parts of a single complex Bose field or, to their independent 
combinations ^ and 

In Eqs. and (jSji, 7c denotes time ordering along a closed 
time-path contour C appearing in the source term integrals 
with = f^dxQ J (Px 1I73I I74I1 . The time path C extends 
from the initial time to to some finite time t > to, and, back, 
from t to to- Note, that in Eq. (|5}, the real-time contour must 
contain the times of interest, i.e., a;o and ya. In practice, this 
is no problem since the largest time, max(2;o, ya), is kept as 
a variable which evolves in the time evolution equations for 
the correlators dicussed below. The second half of C, from t 
to to, ensures the normalization of the generating functional 
Z[0, 0; Pd] = 1, i.e., unity of the trace of the density matrix. 

The generating functional Z[J, K; po] has a functional in- 
tegral representation, which can be found by inserting, to 
the left and right of po, a complete set of eigenstates of 
the Heisenberg field operators at the initial time, $(xo ~ 
to,x)|$±) = $±(x)|$±), with$±(x) = $(to,x): 



Z[J,K-pd]^ j d<^>+d<^-{<^>+\pD{h)\<^- 



x I I?$exp|i[s'[$] + / ^{x)J{x) 
^ ^ ^{x)K[x,y)^{y)\y (6) 



The functional integral J I?$ — J| ■ J sums over all field 
configurations. S denotes the classical action, defined as 



sm = / c{x), 



(7) 



where C is the Lagrangian density. We will provide be- 
low. 

A general initial density matrix can be parametrized as 



($+|pz,(to)|$") =AAcxp{i/c[$]}, 



(8) 



with a normalization A/^ and the functional / expanded in pow- 
ers of the fields: 

/c[4>]=ao + ^— / an{xi,...,Xn)Y[Hx,). (9) 

n=l ' "'2:1. ..a;„ 



Here, the coefficients an{xi, Xn) vanish identically for all 
times different from to (cf., e.g., Ref. 117 Itl ). 

In many practical cases it is sufficient to specify, at time to, 
only the lowest correlation functions. If an initial state is fully 
determined by the mean field 

^,ix) = (10) 

and the connected two-point function or cumulant 

then the initial density matrix, Eq. (|8|i, can be written as a 
Gaussian in the field $, i.e., all with n > 3 vanish identi- 
cally. 

In this case, comparing Eqs. (|8} and (|9} with Eq. (|6jl, one 
finds that the initial-time sources can be absorbed into the 
functional integral by redefining the source fields at times 
2^0 = 2/0 = ^0 according to J{x) — + J{x) — ai{x) and 
K{x,y) — > K{x,y) — a2{x,y). ao yields an irrelevant 
normalization constant. In this way we arrive at the non- 
equilibrium generating functional in the form 

Z[J, K]= jv^ exp S'[$] + ^ ^{x)J{x) 

^{x)K{x,y)^{y)\y (12) 



We emphasize that the use of a Gaussian initial density matrix 
only restricts the "experimental" setup described by the initial 
conditions for correlation functions — higher irreducible cor- 
relations can build up corresponding to a non-Gaussian den- 
sity matrix for times t > to. Non-Gaussian initial density 
matrices pose no principal problems but require taking into 
account additional initial-time source fields. 

The 2PI effective action is obtained below from a double 
Legendre transform of the generating functional M2\ with 
respect to the linear and bilinear source terms, J{x) and 
K{x,y). Accordingly, more complicated initial conditions 
involving higher irreducible correlation functions are most 
efficiently described in terms of nPI effective actions with 
n > 2 1 75]. In the following we will always assume that, 
in order to fully specify the initial state of the many -body sys- 
tem, it is sufficient to provide the initial one- and two-point 
functions 6 and G. 



B. 2PI effective action 

The 2PI effective action is defined as a Leg- 

endre transform of the generating functional of connected 
Greens functions W[J, K], defined by 



Z[J,K] = exp{iW[J,K]}. 



(13) 
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The double Legendre transform reads 



with the proper self-energy 



T[<I,,G]=W[J,K]- J <p{x)J{x) 

(l){x)K{x,y)(j){y) 



xy 



G{x,y)K{x,y). 



xy 



Here, J and K are fixed through the conditions 

5W[J, K] 



6Ji{x) 



(14) 

(15) 
(16) 



with the mean field (t)i{x) and the two-point function 
Gij {x, y) defined in Eqs. (I10> and ([^, respectively. 

The equations of motion for (j)i{x) and Gij {x, y) are given 
by the stationarity requirements 



S(j){x) 
SGix,y) 



0, 



0, 



(17) 
(18) 



which follow directly from the definition J14> for vanishing 
sources J and K. 

Eq. ( I14t shows that the conventional IPI effective action 
r[cj)] is merely r[0, G] for K — 0, i.e., it is equivalent to 
r[(/), G] for that function G for which the stationarity condi- 
tion ( list is fulfilled. Hence, the equation of motion following 
from the condition ( I17> is equivalent to that derived from the 
IPI effective action. On the exact level all effective actions 
are, of course, equivalent. Here the particular reason for us- 
ing the 2PI effective action is that it efficiently allows to de- 
vise systematic approximation schemes suitable for nonequi- 
librium dynamics. 

It is convenient to write the 2PI effective action as 



r[(/.,G] = 5M 



-Tr jlnG" 
2 



const. 



-G,'[cf>]G}+T2[<f>,G] 
(19) 



which contains the contribution from the classical action S, a 
one-loop-type term and a term r2 [0, G] that contains all the 
rest. The trace, the logarithm and the product of Greens func- 
tions in the third term are meant in the functional sense. Gq ^ 
is the inverse of the classical propagator 



^Go}jix,y;(l)) 



5(t)i{x)5(t)j{y)' 



(20) 



for which we will give an explicit expression in the next sec- 
tion. 

Taking the derivative of M9\ in terms of G one observes 
that the second stationarity condition ( I18t is equivalent to the 
exact Dyson-Schwinger equation for the propagator. 



G,Hx,y) = Gg ■ (x, (/i) - Sy(x, y; (/», G), 



(21) 



Y.ij{x,y; 4), G) = 2i 



6Gij{x,y)' 



(22) 



From Eq. (I21> one easily sees that the proper self-energy E is 
IPI. Since the functional differentiation of r2 with respect to 
G corresponds to opening one propagator line in any diagram 
contributing to the expansion of r2 may only contain dia- 
grams which are at least two-particle irreducible. This feature 
constitutes the name of the 2PI effective action. 

We note that the above mentioned, more general nPI effec- 
tive actions can be constructed analogously. These depend on 
correlation functions up to order n and generate diagrammatic 
expansions in which also the corresponding higher vertices are 
self-consistently determined, i.e. dressed |71]. 

Before we apply the effective-action description of a non- 
equilibrium many-body system to derive dynamic equations 
we close this section by introducing the spectral. 



p,j{x,y) ^ i{[<Pi{x),^^{y)])^, 
and statistical components, 

F,,{x,y) = h{^,{x),^,iy)}),., 



(23) 



(24) 



of the two-point function G, where the subscript c means that, 
for 7^ 0, disconnected parts are substracted. In Eq. i24\ . 
{, } denotes the anticommutator The decomposition identity 
for the two-point function G follows from Eqs. (II 1> . ( I23> . and 
(ES as 

Gtj{x,y) = F^j{x,y) - -py (a;, y) signc(a;o - yo), (25) 

where sign^(a;o — j/o) is the sign function along the path C 
which evaluates to 1 (—1) if xq (yo) is prior to yo (xq) along 
C. 

While the spectral function encodes the spectrum of the the- 
ory, the statistical propagator gives information about occupa- 
tion numbers. Loosely speaking, the decomposition makes 
explicit what states are available and how often they are occu- 
pied. Note that, in thermal equilibrium, the spectral and statis- 
tical functions would be related by the fluctuation-dissipation 
relation, which, for a non-condensed homogeneous system in 
energy-momentum space, reads 



'I I ^+n{uj,T) ) p 



(cq) 



(^,P), (26) 



with the Bose-Einstein distribution function n{uo,T) — 
(^g('^-p-)/kBT _ Relation (l26t does no longer apply in 

a system far from equilibrium, such that the dynamical evo- 
lution of F and p is given by separate (coupled) equations of 
motion. 



C. Dynamic equations 

In this section, we derive dynamic equations for the mean 
field (f) and the two-point cumulant G from the stationarity 
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conditions, Eqs. ( I17> and (I18> . We consider explicitly the the- 
ory defined by the Hamiltonian Q and choose the representa- 
tion 



$i(x) = \/2Re*(a 



$2(x) = V2Im*(x) (27) 



of the quantum field in terms of its real and imaginary parts. 
The corresponding field operators obey the bosonic commuta- 
tion relation [$i(x, t), $2(y, t)] — i<5(x — y) while all other 
equal-time commutators vanish. This representation will be 
particularly convenient when discussing the non-perturbative 
expansion in Chapter^] The classical action is given 

by 



1 



M^)tD-^\x-y)<f,{y) + Sint [$] , (28) 



with the inverse free propagator (Greens function) 

"i^ij^i^ -y) ^Sc{x- y) {-ia2,tjdy„ - iJiB(x)(5„) , 

(29) 



containing the Pauli matrix 

0-2 = 



-i 

1 



(30) 



Here, Sc{x — y) = Sc{xq — yo)S'^^'>(x — y) denotes the four- 
dimensional Dirac distribution on the closed time path. The 
interaction part 5int [^] may be expressed, for the Hamiltonian 
(□, with V{x -y) = F(|x - y|, xo)^^,^;,, as 



^i„t[*] 



t>,{x)^,{y)V{x - y)^,{y)^,{x). (31) 



xy 



where it is summed over double indices. We use Eqs. ( 128] 
to derive the classical inverse propagator, cf. Eq. \2Q\ : 



5(t)i{x)5(t)j{y) 



^DiMx -y)~ (j)i{x)V{x - y) (f)j{y) 



(t>k{z)V{x - z) (I)k{z)5ij5c(x - y), (32) 



Now, everything is prepared to derive the dynamic equations 
by imposing the stationarity conditions, Eqs. ( I17> and ( I18t . 
These conditions yield, together with Eqs. il9\ . ( I28> . ( I31l l. 
and(|32}: 



= ^[Goj,:(y,^,0)-G7/(2/,x)] 
+ S^2[c^,G]/SG^Jix,y). 



(34) 



From Eq. i33\ . the equation of motion for the mean field is 
obtained as 

[ - i(T2,tjdxa - i?lB(x)%] 

V{x -y) ^ (^fe {y)<l>k (y) + Gkk {y, y)) {x) 



+ Gtj{x,y)(l)j{y) 
ST2[(j),G]/SMx). 



(35) 



The dynamic equation for the Greens function G follows from 
Eq. ( I34> by convolution with G and using the definition (I22> 
of the self-energy E: 



Go}ki^^z;(f>)Gkjiz,y) 

6c{x-y)Sij+ / T.ik{x,z;(j),G)Gkj{z,y). (36) 



Once the 2PI part T2[4>, G] of the effective action is known, 
Eqs. i35\ and i36\ . representing a closed system of dynamic 
equations for (j) and G, may be solved to determine, for given 
initial correlation functions 0(to) and G{to, to), the complete 
dynamical evolution of the system. It is, of course, not possi- 
ble to solve the equations exactly. In the following we discuss 
approximation schemes for r2 and use these to obtain approx- 
imative dynamic equations for and G. 



D. Number conservation 

A major advantage of the effective action approach is that 
it automatically provides us with dynamic equations which 
are particle number conserving. This is a consequence of the 
Noether theorem in conjunction with the invariance of the the- 
ory under orthogonal transformations and can be seen as fol- 
lows lloill : The stationarity conditions (I17> and (I18> for the 
mean field and the propagator can be combined to the equa- 
tion 



S(j)j{x) 



ST[(t),G] 
SGkj (y, X 



■Gki {y,x) 



0, 

(37) 



o = sr[cb,G]/SMx) 

= f \tD^'ix-y)4>,ix) 

- V{x ~ y)^(^(pk{y)<Pk{y) + Gkk{y,y)^4>t{x) 

-V{x-y) G.i-i{x,y)(t)j{y) 
+ 6r2[^,G]/6Mx), (33) Here, 

= Jr[0, G]/6G,j{x,y) 



with the elements a2.ij of the Pauli matrix i30\ . From the 
specific expression ( I19> for the 2PI effective action follows 
that Eq. i37\ is equivalent to the relation 



dxon{x) - Vj(.T) 



= —2ia2.'. 



4>t{x) 



,G] 



S(j)j{x) 



ST-,. 



,G] 



SGkj{y,x) 



Gki{y,x) 



n(x) = (l)i(x)(l)i{x) +Gii{x,x), 
i{x) — [02(a;)V(/)i(a;) -(/)i(a;)V<?!)2(a;) 



(38) 



(39) 



(rc($2(x)V$i(x) - $l(x)V$2(x)))c] (40) 



are the total number and current densities, respectively. 
Clearly, particle number is conserved locally if n and j obey a 
continuity equation, i.e., if the right hand side of Eq. (I38> van- 
ishes identically. We consider the specific structure of these 
terms: The interaction part of the 2PI effective action occur- 
ing therein is defined as 



ri„t[0,G] = r[0,G] 



Mx)iD^j^{x - y)(f>j{y) 



^Tr [D-'G] 



(41) 



The 2PI effective action is, like the underlying quantum ac- 
tion S[^], Eqs. j28> and i31\ . a singlet under 0(2) rotations. 
It is parametrized by the classical fields ipi and Gij, where 
the number of 0-fields has to be even in order to construct an 
0(2)-singlet. From the fields (pi alone one can construct only 
one independent invariant under 0(2) rotations, which can be 
taken as tr{(p(j>) = 0^ = (picpi. All functions of 4> and G, which 
are singlets under 0(2), can be built from the irreducible, i.e., 
in field-index space not factorizable, invariants 1 3 3, . .34] 



tr(G"), 



and 



ti-(</)(/.G"), 



(42) 



with n — 1,2, .... As before, the trace tr(-) only applies to 
the field-component indices while there is no integration over 
space-time, e.g., tr(G^) = Gij{x,y)Gjk{y,z)Gki{z,x). 

For contributions to Tint which contain only 0^ or tr(G"), 
the terms in square brackets in Eq. (I38> either vanish sepa- 
rately or are symmetric under the exchange of i and j. More- 
over, if a term contains an invariant of the form tr((/)0G"), 
as, e.g., the contributions remaining in Tint from Tr{GQ ^G}, 
the combination of the terms in square brackets in Eq. (I38> 
is symmetric under transposition in field index space. Hence, 
the total number density is conserved locally as a consequence 
of the 0(2) symmetry of the theory, and, more importantly, 
this is true for any set of approximative dynamic equations 
derived from a truncated but still 0(2)-symmetric effective 
action. Note, finally, that only terms in the action which con- 
tain mixed invariants tr(00G" ) induce exchange of particles 
between the condensate and the non-condensed fraction of the 
gas. 



III. DIAGRAMMATIC REPRESENTATION 
A. Loop expansion of r2 [0, G] 

In this section we briefly discuss the diagrammatic meth- 
ods that we will use later on to generate the equations of mo- 
tion for (j) and G. For a more comprehensive review see, e.g., 
Ref. fnll and references cited therein. 

The 2PI effective action, Eq. ( I19> . can be expanded into 
a series of closed loop diagrams, i.e., diagrams without ex- 
ternal lines. In Fig. [2 the diagrammatic representation of 
r[<j), G] - T2[<j), G] is shown, cf. Eq. (|19}. Fig. 13 shows the 
lowest order diagrams contributing to r2[(j>, G]. Internal lines 
in these diagrams represent the full two-point Greens func- 
tion G, vertices the bare interaction V. In addition, one needs 



r "-P[(^,G] = + ]^wwv^ + i Tr In G 




FIG. 1: (Color online) Diagrammatic representation of the one-loop 
part r"°°P[<j!-, G] = r[<;/), G] - r2[0, G] of the 2PI effective action, 
cf. Eq. <19t . without the irrelevant constant. The inverse free prop- 
agator D~^(x,y) is represented as a (green) filled circle, external 
fields (t>{x) as (red) crosses, the bare potential V{x — y) as a wiggly 
line, and the full propagator G{x, y ) as a (blue) solid line - colors are 
seen in the online version. At each 'vertex' it is summed over double 
field indices and integrated over the respective space-time variable 
X. The first two diagrams are the classical action the last three 
diagrams the term |Tr{Go ^G}, cf. Eq. 02l . 



r231°°P[(^,G] 



-I- 



FIG. 2: (Color online) Diagrammatic representation of the two- and 
three-loop diagrams contributing to the 2PI part r2 [<j}, G] of the 2PI 
effective action, cf. Eq. <19> . To simplify the diagrams, the bare ver- 
tices V{x — y) are drawn as black dots. Each such vertex is under- 
stood to represent a sum of the topologically different terms shown 
in Fig.|3| At each vertex, it is summed over double field indices and 
double space-time variables according to the respective diagram in 
Fig.El 



to distinguish topologically different diagrams. This becomes 
obvious when considering an interaction V = V{x — y, xq) 
which is non-local in its spatial arguments and its diagonal 
coupling to the field components 1 and 2 as exhibited by 
Eq. ( 13 U . An example of this are the diagrams shown in Fig.|3 
which combine to give the "double-bubble" diagram, i.e., the 
first graph in the loop expansion of r2 in Fig.|2] 

Counting the number of loops, r[<j), G]~T2[(j), G] is of lead- 
ing order, and, as discussed in Appendix IaI this part of the 




V(x-y) = = I -I- ^(\AA/\^ 





FIG. 3: The representation of the bare vertex in terms of a black dot 
stands for a sum of the three topologically different connections of 
the four 'comers'. The black lines do not represent propagators and 
are only drawn in order to illustrate the different possible connections 
of propagators and/or external fields at the vertices in Figs.|4|and|5| 
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r2"pB[(i),G] = 



FIG. 4: (Color online) Diagrammatic representation of the diagrams 
contributing, in the Hartree-Fock;(-Bogoliubov) approximation, to 
the 2PI part r2 [(t>, G] of the 2PI effective action. At each vertex, 
it is summed over double field indices which are indicated, and in- 
tegrated over double space-time variables. For an explanation of the 
symbols see the caption of Fig. Q 



effective action equals the IPI effective action to one-loop or- 
der. We have seen in Section HlBl that r2 can only contain 2PI 
diagrams. For this reason it is at least of two-loop order, the 
leading order comprising the double-bubble and the "setting- 
sun" diagrams in Fig.|2] 

In the following section we will consider the most sim- 
ple approximation of r2 which only takes into account the 
double-bubble diagram. The dynamic equations resulting in 
this approximation constitute the well known time-dependent 
Hartree-Fock (vanishing mean field, (j) ~ 0) or Hartree-Fock- 
Bogoliubov {(j) 7^ 0) equations which describe coupling be- 
tween the condensate and the thermal fractions of the gas and 
account for the formation of pair correlations through binary 
interactions. 



B. Hartree-Fock-Bogoliubov dynamics 

In this section we derive the coupled Hartree-Fock- 
Bogoliubov equations of motion for the mean field and the 
correlation function G by evaluating Eqs. i35\ and i36\ for the 
leading-order approximation to the 2PI part of the effective ac- 
tion, the double-bubble diagram in Fig.|2]which is shown, for 
a non-local coupling V{x — y), in Fig.0] 

For a better comparison with standard results, we shall, in 
this section, use the representation where the field $ and its 
complex conjugate are considered as independent quantities. 
To this end we write the effective action in terms of a two- 
component field with 



$i(x) = *(a;), $2(2;) = **(a;), 



(43) 



Note that the representation in terms of the components de- 
fined in ( I43> is obtained from the 'cartesian' representation 
(I27> by the U{2) unitary transformation 



$1 
$2 



J_ / 1 i 
^/2 U 



$1 
$2 



(44) 



Hence, introducing an upper-index notation according to 

G^^^^^Cj, etc., (45) 
with i,j £ {1, 2}, the double-bubble contribution to r2 can 



be written in the form 



rf«[<^,G]=~- V{x-y) G:{x,x)G^\y,y) 



xy 



+ 2G^\x,y)G;{y,x) 



(46) 



where it is again summed over double indices. The statistical 
factor 2 in the second term in brackets accounts for the topo- 
logically equal diagrams, i.e., for the second and third diagram 
in Fig. El 

The proper self energy corresponding to the F^^^ in 



Eq. (I46> is obtained using Eq. (I22> : 

ST2\<P,G] 



E/(x,2/) = 2j 



5G,'[x,y) 
i 
2 



Viz-x)G,'{z,z)S/5c{x-y) 



+ 2V{x-y)G;{y,x)y 
By making use of the symmetry relations 

F^\x,y)^F^,{y,x)^[F/{y,x)]*, 
P'{x,y) = -p\{y,x) = [-pHy,x)]* 



(47) 



(48) 
(49) 



the dynamic equation \2>5\ for the mean field may be reex- 
pressed in terms of (j) — (fii = (fp' and F as 



= / y(x-y,t) 



(y,i)0(x,i)+F„,(y,x;O0*(y,t) 



+ F„(y, x; t) <^(y, t) + F„ {y, y; t) </)(x, t) 



(50) 



In deriving Eq. i50\ we have integrated yo over the path C, 
leaving only the space integral ~ J d^y, and used the lo- 
cality of V{x — y) — V{x — y, t)dc{xo — yo) in its time ar- 
gument t = xq. By this integration, the spectral functions 
are evaluated at equal times and disappear from the equa- 
tion, either because of a vanishing equal-time commutator, or 
through cancellation between different terms. The statistical 
two-point functions are also evaluated at equal times, and we 
denote the different matrix elements occuring as 

F„(x,y;t) = F,'{x,y)\^^^^^^ = i({*t(x, t), *(y, t)})„ 

(51) 



^^„(x,y;t) = F^\x,y)\ 



-({*(x,i),*(x,t)}). 



(52) 



Eq. iSOi is the well-known non-linear dynamic equation 
for the condensate mean field, where the contribution from 
the three-point connected Greens function or cumulant 
(#t(y)4((y)#(x))^ has been neglected, cf., e.g., Ref. |68i|. 
The equation, however, involves the symmetrized two-point 
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cumulants instead of the normal ordered ones, the time- 
dependent non-condensate and anomalous density matrices, 
n and rh, respectively: 

n(x,y;t) - (*t(y,t)^(x,<)),, (53) 

m(x,y;t) = (*(y,i)*(x,t))e. (54) 

The reason for this is that the functional integral Z, as defined 
in Eqs. il2\ . ( I28> , is equivalent to a Hilbert-space formula- 
tion with a Weyl-ordered Hamiltonian operator H, i.e., an 
H which is invariant under permutations of non-commuting 
field operators. Hence, the effective-action approach yields 
dynamic equations involving symmetrized correlation func- 
tions Fn and F,„. In order to rewrite Eq. ( 15 Ot in terms of 
normal ordered cumulants one can use the commutation rela- 
tions which give F„(x, y; t) — n(x, y; t) + S{x — y)/2 and 
F,„(x,y;t) = to(x, y;t). This introduces, in Eq. (l5Qt . infi- 
nite terms involving the combination V{0, t)+ V(y, t)d{0), 
which also appear if the mean-field equation is derived by 
means of the Schrodinger equation for the density matrix in 
Eq. (|3|l, with a Weyl-ordered Hamiltonian. These contribu- 
tions disappear when shifting the zero of the energy scale 
through normal ordering of the Hamiltonian. 

The mean-field equation i50i shows that, in order to ob- 
tain a closed set of dynamic equations, the equations of mo- 
tion for Fn{x,y) and Fm{x,y) are required only for equal 
times, xq = yo. The equation for Fn{x, y) is obtained from 
the equivalent of Eq. ( I36t in the representation i43l . by sub- 
stracting the expression for the time derivative dygG2{y, x) 
from that for dx„Gi{x, y). Taking the equal-time limit L921 
2^0 = 2/0 ^ 0, the resulting equation reads 

[ihdt ~ i?iB(x) + HiB(y)] F„(x,y;t) 



y(x-z,i) 

(^(x, t)cj){z, t) + F„, (x, z; t)) F,;(z, y; t) 
+ r (z, t) ((/.(z, t)i^„(x, y; t) + 0(x, t)F„(z, y; t) 
+ F„(z,z;i)i^„(x,y;<) + F„(x,z;i)F„(z,y;t) } 



- {x y}* 



(55) 



The last term denotes the complex conjugate of the first term 
in curly brackets, with x and y interchanged. Similarly, the 
time evolution of i^m(x, y\t) is derived by adding the 

\iovisfordxoGi{x,y) andfor9yo(G'2^(y,2:))*, and taking the 
same equal time limit as above: 

[ihdt ~ i/iB(x) - i/iB(y)] Fm(x,y;t) 

y(x-z,t) 

X (^0(x,t)0(z,i) + F™(x,z;t)ji^„(y,z;<) 

+ 0*(z, t) ((/)(z, <)F,„(x, y; t) + 0(x, t)F,„(z, y; i)) 

+ F„(z,z;i)F™(x,y;t) + K(x,z;i)Kn(z,y;t) | 
{x ^ y}. (56) 



Reexpressing Eqs. J55> and \56\ in terms of normal ordered 
cumulants n and rh |93], one recovers, together with the 
non-linear mean-field equation ( 15 Ot . the familiar closed set of 
dynamical Hartree-Fock-Bogoliubov equations describing the 
exchange of atoms between the condensate and thermal frac- 
tions of the gas; cf., e.g., Ref. i68ll . These equations have 
been extensively used to describe the short-time, far-from- 
equilibrium dynamics of strongly interacting Bose-Einstein 
condensates. Applications include the production of corre- 
lated atom pairs in four- wave mixing of condensates |68|, the 
atom-molecule coherence in a Ramsey-type interferometer 
1 69, 78], the condensate and molecular dynamics in Feshbach 
crossing experiments 179, .80. .8 LI . as well as the many-body 
dynamics near a photoassociative resonance 1 82, 83]. 

Note that the "diagonal" part of the equation of motion ( I56l l 
for Fjn, i e-, the term in square brackets on the left-hand side, 
receives an additional term y (x — y, t) fri{x, y; t) when the 
equation is reexpressed in terms of the density matrices h and 
rh: 

[ihdt - HiB (x) - HiB{y) - V{x ~ y, t)] m(x, y; t) 
= y(x-y,t)0(x,t)0(y,t) + ... (57) 

Hence, the diagonal part of the dynamic equation for rh con- 
tains the two-body Hamiltonian H2B{x,y) = Hib{x) + 
Hmiy) + V{x — y) instead of the free Hamiltonian without 
V . It was shown in Ref. |68] that this is crucial when deriv- 
ing the Gross-Pitaevskii equation (GPE) from a many-body 
Hamiltonian with a non-local coupling V{x — y), in the limit 
where the gas is weakly interacting, i.e., dilute, and the inter- 
action strength g is constant in time. In the GPE, this coupling 
constant g multiplies the term non-linear in </) and is defined 
in terms of the s-wave scattering length a as g ~ Airh^a/m. 
The reason for this is that g is equal to the zero-energy limit 
g — (277/1)"^ limp_^o(p|72B(p^/"J + *0)|p) of the quantum- 
mechanical two-body transition (r)-matrix 



(p|r2B(i?)|p') - {v\v{i + G2b{e)v)\v') 



(58) 



Without the term proportional V in Eq. \51l . the full two- 
body Greens function G2B = [E — H2b)^^ in the T-matrix 
( I58l l would be replaced by the free Greens function G'2b.o = 
{E — H2B + y)^^, see Ref. [68]. As a consequence, the cou- 
pling g involves only the second-order Born approximation 
T|bI^°™(£;) = V{1 + G2B,a{E)V) of the T-matrix, which, 
in general, does not give the s-wave scattering length in the 
zero-energy limit. This also shows that the non-local bare cou- 
plings V{r) at different interatomic distances r do not form a 
suitable set of small parameters for a perturbative expansion. 

In summary, taking into account the double-bubble dia- 
gram as the single contribution to yields the closed set of 
Hartree-Fock-Bogoliubov many-body equations for normal- 
ordered equal-time one- and two-point correlation functions, 
i.e., the mean field and the normal and anomalous density 
matrices. Moreover, starting with a non-local (bare) cou- 
pling V{x — y), the full GPE requires contributions from the 
double-bubble diagram, i.e., beyond one-loop order On the 
other side, when starting with a local coupling V{x — y) = 
g6{x — y), the underlying assumption is that the Hamiltonian 
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defines an effective theory which is valid only at low colli- 
sion energies. In such a theory, the GPE directly results as the 
"tree-level" approximation, where all contributions from two- 
and higher 71-point functions are neglected. 

We finally remark that the spectral functions do not 
appear in the HFB equations for the equal-time correlators. 
Hence, the HFB equations remain invariant for vanishing 
commutators of the Rose fields and can therefore be regarded 
as fully classical. 



IV. NON-PERTURBATIVE EXPANSION 

The dynamics of non-dilute ultracold Bose gases can not, in 
general, be described using a loop or coupling expansion. In 
three spatial dimensions, this is the case, for instance, when 
the mean interparticle distance is considerably smaller than 
the scattering length, r] = na^ ^ 1. This means that the 
approximation needs to be based on an expansion parameter 
different from the interaction strength. 

An expansion in powers of the inverse 1/M of the number 
of field components provides a controlled expansion parame- 
ter which is not based on weak couplings. It can be applied to 
describe physics characterized by large fluctuations, such as 
encountered near second-order phase transitions |60|, or for 
extreme non-equilibrium phenomena such as parametric reso- 
nance \i5>\ or spinodal decomposition dynamics |38]. For the 
latter cases a 2PI coupling or loop expansion is not applicable. 

In the following sections we develop such a scheme for a 
single-component Bose-gas, on the basis of an expansion in 
inverse powers of the number A/^ of field components. For de- 
tailed studies of such an expansion in the context of relativistic 
quantum field theory confer Refs. L33, .34.1 . The method can 
be applied to bosonic or fermionic theories alike if a suitable 
field number parameter is available, and we exemplify it here 
for the case of a scalar 0(A/')-symmetric theory. Again, we 
point out that the apparently rapid convergence of the 2PI ex- 
pansion even for small values of M ~ 2 is crucial for our 
approach |60, 61]. 

We will first provide the relevant terms in a expansion 
of the 2PI effective action to next-to-leading order and use 
these to derive the dynamic equations. 



A. 1 /M expansion of the 2PI effective action 

In the following we consider a systematic non-perturbative 
approximation scheme for the 2PI effective action. It classi- 
fies its contributions according to their scaling with powers of 
1/N, where M denotes the number of field components. We 
will assume that no further U{1) symmetry characterizes the 
theory, i.e., we choose all M field components to be real. 
For A/^ = 2 the 0(2)-symmetric theory for real fields which 
we again denote as $1 and $2, is equivalent to the U{1)- 
symmetric theory considered in Section Ull] for the complex 
fields- = ($1 +i$2)/V2. 

The 2PI 1/M approximation scheme has been derived for 
a relativistic field theory with local interactions in Refs. |33l 



13411 . Before we proceed, we would like to point out that the 
2PI 1/M scheme is different from the standard IPI re- 
summation which involves the free instead of the full 2-point 
Greens function. Here we consider a non-relativistic field the- 
ory with the classical action given by Eq. (I28> . Recall that 
connected n-point Greens functions with n > 2 are equivalent 
to ordinary n-point functions for the fluctuation field ^i{x) 
obtained by substracting the mean field from the field opera- 
tor, ^i{x) = (t>i{x) + ^i{x). Hence, the diagrammatic expan- 
sion of the 2PI effective action part r2 can be derived by use 
of the non-local effective interaction {y{x ~ y) = V{y — x)) 

+ $,(x)l>,(y)l>,(2/)$,(a:)l, (59) 



which is obtained from Eq. (13 1> by shifting $i (x) (pi (x) + 
^i{x) and collecting all terms cubic and quartic in the fluc- 
tuating field ^i{x). One observes that, in the presence of a 
non-vanishing macroscopic field 4)i{x), there exists an effec- 
tive three-vertex in addition to the four-vertex. 

The 1 /M classification scheme is based on invariants un- 
der 0{M) rotations which parametrize the 2PI diagrams con- 
tributing to r[(p, G]. The 2PI effective action is, like the under- 
lying quantum action a singlet under 0{M) rotations. 
As mentioned already in Section HTdI only the 0{M) invari- 
ant Xicfxp = (f)^ — (f)i(t)i ^ M can be constructed from the 
mean field alone. All other functions of (j) and G, which are 
singlets under 0{M), can be built from the irreducible invari- 
ants given in Eq. ( I42> above. From these, for a given M, only 
the invariants with n < M aie linearly independent — there 
can not be more invariants than fields |94|. In particular, for 
the next-to-leading order approximation one finds that only 
invariants with n < 2 appear, which makes it very suitable for 
practical calculations. 

The factors of A/^ in a given graph contributing to r[0, G] 
have two origins: Each irreducible invariant is taken to scale 
proportional to M since it contains exactly one trace over the 
field indices, while each vertex provides a factor of 1/M, 
cf. Eq. ( I59> . The factors of M arising from the irreducible 
invariants (142 > correspond to the number of closed lines fol- 
lowing the field indices in a diagram, plus the number of lines 
connecting two insertions of the classical field 4>. 

The expression (I19> for the 2PI effective action contains, 
besides the classical action, the one-loop contribution propor- 
tional to Tr (In G-^ + Go^[<j)]G), and a T2 [<j), G] which con- 
tributes if higher loops are taken into account. The one-loop 
term contains both leading order (LO) and next-to-leading- 
order (NLO) contributions in an expansion in powers of 1 /M. 
The logarithmic term corresponds, in absence of other terms, 
to the free-field effective action and scales proportional to 
M. To separate the LO and NLO contributions at the one- 
loop level consider the second term Tr {GQ^[(t>]G). From the 
form of the classical propagator, cf. Eqs. (I29> . ( I32> . one ob- 
serves that it can be decomposed into a term proportional to 
tr(G) ~ M and terms - {V/M)[ti{<f)(t))ti{G) + 2ti{(f><f)G)]. 
This can be seen as the sum of two "2PI one-loop graphs" with 



field insertion 



and 



J, respectively, as shown in 
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(a) 

r2i-o[(|),G] = 

(b) 

r,'^LO[(|, Gl : 




FIG. 5: (Color online) Diagrammatic representation of the diagrams 
contributing, in leading order (LO) and next-to-leading order (NLO) 
of the 1 /A/'-expansion, to the 2P1 part r2 [(p, G] of the 2P1 effective 
action, cf. Eqs. I60i and j6H . At each vertex, it is summed over 
double field indices and integrated over double space-time variables. 
For an explanation of the symbols see the caption of Fig.Q 



Fig H Counting the factors of J\f coming from the traces and 
the prefactor, one finds that only the term ~ tr((/)(/)) tr(G) con- 
tributes at LO, i.e., it scales proportional to JV, while the term 

~ tT{(t)(j)G) is of NLO. 

The LO contribution to T2[(f>, G] consists of only one two- 
loop graph, which is ^-independent and corresponds to the 
first contribution to Vf^^ in Eq. (|46} (cf. Fig.|4j: 



r2°[G] = / Gu{x,x)V{x - y)G,,{y,y) . (60) 



xy 



It is shown diagrammatically in Fig. [5^. In NLO (Fig. [Sji) 
there is an infinite series of contributions which can be 
summed up analytically l33ll34ll : 



r^^"['^,G]=.-Trlii(i?[G]) 



I{x, z; G)V{z - y) (j),{x) G^^j {x, y) (j)j{y) 



(61) 



FIG. 6: (Color online) Recursive definition of the resummed local 
interaction function A(3::, y; G), cf. Eqs. {65\ and <64L in NLO of 
the 1 /J\f expansion. The thick wiggly line represents A, while G is, 
again, depicted as a thiclc (blue) line, and V{x — y) a&2L thin wiggly 
line. Colors are seen in the online version. 



with 



Tr HB{G)] 



xy 



X ( i-^ / V{y - w)Gkiiw,x)Gki{w,x) 



(63) 



Iix,y;G) ^ jj: I V{x - z)G^j{z,y)Gt.j{z,y) 



ijr / I{x,w;G)V{w ~ z)Gij{z,y)Gij{z,y) 



AAA 



A(a;, z; G)Gij (z, y)Gi.j (z, y), 



(64) 



A{x, y; G)=V{x-y)-i I{x, z; G) V{z - y), (65) 



The diagrammatic representation of A is shown in Fig.|6] The 
first term on the rhs of Eq. i63\ corresponds to the two-loop 
graph with the index structure exhibiting one trace such that 
the contribution scales as Xi{G'^)/JV ~ Af^. This graph is the 
first in the expansion of F^'^'-' shown in Fig.lSji. One observes 
that each additional contribution scales as well proportional to 
[tr(G2)/7V]" ^ TV" for all n > 2. Thus all terms contribute 
at the same order. 

The functions I{x, y; G) and the inverse of B{x, y; G) are 
related by iS 13 



B-\x,y;G) = Six - y) - iI{x,y;G), 



(66) 



which follows from convoluting J62t with B ^ and using 
Eq. (I64> . We note that B and / do not depend on (j>, and 



Here, 



B{x,y;G)=6{x-y) 



(67) 



+ I Vix- z)G,j{z,y)G,,iz,y), (62) 



r2[(/.,G] = r^"[G] + r^^"[0,G] 



is only quadratic in (p at NLO. It was shown in Ref. 116 IB 
that invariants containing more than two three-vertices, i.e., 
field insertions, are two-particle reducible and therefore can 
not contribute to Furthermore, it was shown that graphs 
containing other invariants, e.g. tr(G"^) are of higher order 
~ 0{1/M). 
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B. Dynamic equations 



Here, the function P is defined as 



In this section we employ the stationarity conditions ( I17> 
and ( I18> to derive the equations of motion for the mean field 
(f) and the two-point function G to NLO in the 1 /Af expan- 
sion of the 2P1 effective action, Eq. ( I67> . In the following, 
we restrict ourselves again to a single complex Bose-field, 
Af = 2, which, in order to simplify the equations, we describe 
by its rescaled real and imaginary parts <& = ($1 + 1^2) / V^, 
Eq. int . In this representation all fields occuring in the equa- 
tions will be real-valued. We can take over all equations 
from earlier sections by lowering all indices and identifying 
(j) = (pi, (f)* = (j)2- The equation of motion for 4)i{x) can then 
be read off Eqs. (I35> and i6l\ : 

[ - i<y2Sjdxo - -fflB(x)] (t)i{x) 

= \jj{x~y)[(<l>{yf + Gkk {y,y)) c|>^ {x) 

+ {Gij {x,y) + Gji {y, x)j (y) 

I{x, z; G)V{z - y) G,j{x, y)<i>, (y), (68) 

y,z 

where we have used Gij {x, y) — Gji{y, x) and the represen- 
tation of the inverse free propagator given in Eq. i29\ . 

To derive the dynamic equation ( I36l l for Gy in NLO of the 
l/Af expansion, we need 



dl{u, v) 

SG^j{x,y) 



, SB{w,z) ^ 



= - / B-\u,w) V{w-x)B-\y,v) 

+ V{w~y)B-\x,v)\G,j{x,y), (69) 



where we have used Eqs. (I62> and (I66> and suppressed the 
functional argument G in i? and /. We use Eqs. ( I60t . ( I61> and 
\69\ to derive the self energy Ej^ (x, y), Eq. (I22t . in NLO. 

We then express the dynamic equations for (pi and Gij in 
terms of the generalized "mass" matrix 

Mij{x,y](j),G) ^ 5ij5c{x - y) Hib{x) 
I 



(70) 



+ V(.T-z)(0fc(z)(/)fe(z) + Gfefe(z,z)) 

+ Vix~y) (^Mx)<p3 (y) + G,j (.T, y)) 
and the non-local self-energy llo^l 

T.ij{x,y](j),G) = -i[K{x,y;G) - V{x-y)] 

X (l)^{x)(j)j{y) + Gij{x,y) 
- P{x,y\(j),G) Gi-i{x,y). 



Fix, y; 0, G) = / A(x, v; G) H{v, w; 0, G) A{w, y; G), 

J VW 

(72) 

with 

H{x,y;<j),G) ^ (j)iix)Gij{x,y)(l)j{y), (73) 

where, as usual, it is summed over double field indices. 
The dynamic equations then assume the compact form: 



- i<J2,ijdxo - Sij^ I V{x - z)(j>k{z)(l)k{z) (pjix) 



Mij(x,y\(j) = 0,G) +iEij(a;,y;0 = 0; G) (j)j{y), 



(74) 



■i(^2MdxoGk]{x,y) ~ iSc{x ~ y)Sij 

Ahk{x,z;(j),G) +iY.ik{x,z]4>,G) Gkj{z,y). 



(75) 



The dynamic equation (I75> for the propagator, multiplied 
from the left by the free classical propagator D, can be di- 
agrammatically represented as the Dyson series shown in 
Fig. El 

Note that one recovers, from Eqs. ( I74t and ( I75t . the dy- 
namic equations in two-loop approximation of the 2P1 effec- 
tive action by setting A(a:;, y) = V{x — y) in Eq. ( 17 1> . and 
the Hartree-Fock-Bogoliubov equations by additionally set- 
ting P = 0, i.e., by setting E^ = 0. 

To obtain numerically tractable equations of motion we fi- 
nally express all functions in terms of their spectral and statis- 
tical components defined in Eq. (125 > and, for S^j, /, P, and 
H accordingly: 



^^j{x,y) = ^fj{x,y) - -TsP^j{x,y)signc{xo - yo), etc. 

(76) 



Using, furthermore. 



Pij{x,y)dxoSigncixQ - ya)/2 = Pijix,y)Sc{xQ - yo) 

= -io-2,ijSc{x - y) (77) 

(71) we obtain the dynamic equations 
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i<J2,ijdx„ - S^j^ J V{x - z)(pk{z)(l)k{z) (j)j{x)- J Mij{x,y;<j> = 0,F)(j)^{y) = J dy^^{x,y;(j> = Q;G) (j)j{y), 

(78) 



icr2akdxoFkj{x,y) - / Mikix, z;4}, F) Fkj{z,y) 

J z 

icr2,ikdxoPkj{x,y) - / Mik{x, z; (f>, F) pkj{z,y) -- 



dzT.1^{x,z](j},G)Fkj{z,y)- I dzT.^f,{x,z](j),G)pkj{z,y), 
Jo 

(79) 



dzY.1^{x,z;4>,G)pkj{z,y), 



(80) 



Vo 



r 



where dx 

M,,{x,y-(j), F), 



dxo J d^x. The real functions 
{x,y;(f>,G), and T.'lj{x,y;(j),G) are all 
regular in xq and are explicitly provided, in terms of the fields 
(j>i{x), Fij {x, y), and pij {x, y), in AppendixiBl 

Note that, as discussed in Section lTTOl the effective action, 
to any order of the 1/J\f expansion, is manifestly invariant 
under 0{M) rotations, i.e., for A/^ — 2, the dynamic equations 
locally conserve, to any order in the expansion, the particle 
number density. 



NON-EQUILIBRIUM DYNAMICS IN ONE SPATIAL 
DIMENSION 



In the following we will apply the dynamic equations de- 
rived from the 2P1 effective action, to next-to-leading order 
(NLO) in the \/M expansion, to the case of an ultracold 
atomic Bose gas in one spatial dimension. A full dynami- 
cal theory of the formation of a Bose-Einstein condensate out 
of a gas with non-equilibrium distribution of single-particle 
energies, as occurring in everyday experiments, is still an 
open and intriguing problem in the field of atomic matter- 
wave physics |31]. As a first step towards a more complete 
description of the dynamics of crossing the Bose-Einstein 
phase transition we study the time evolution of a uniform 



one-dimensional Bose gas which is, initially, characterized by 
a non-equilibrium distribution of particle momenta. In this 
chapter, we present results for the time evolution of such a 
system as obtained from numerically solving the dynamical 
equations introduced above. 

In one spatial dimension, there is no phase transition to a 
Bose-Einstein condensate, i.e., to a state where only the low- 
est energy mode is macroscopically populated. Consequently, 
there is no condensate mean field. Nevertheless, at very low 
temperatures, equilibrium states with a macroscopic occupa- 
tion of the lowest energy levels are possible, a phenomenon 
which manifests itself in approximate off-diagonal long-range 
order (ODLRO), i.e., in a first-order phase coherence which 
extends over the entire gas sample |7] but can, under specific 
experimental conditions, be broken more easily than in three 
dimensions. On general grounds one expects that, starting 
with a non-equilibrium momentum distribution, the system is 
driven towards a state described by a Bose-Einstein equilib- 
rium distribution with a suitable dispersion relation. 

For vanishing mean field, the NLO dynamical equations for 
the statistical and spectral correlation functions F and p, re- 
spectively, are given in Eqs. (|29j and (|80|l, with </) ee 0, and 
space-time variables x — {xq,xi), etc. We have chosen a 
local coupling 



+ — + 




V{x -y)= gm5{x - y), 



(81) 



where the coupling constant in one dimension, giD, is related 
to the dimensionless 7 parameter and the total line-density of 
atoms ni by 0] 



FIG. 7: (Color online) Diagrammatic representation of the 
Schwinger-Dyson equation for the full propagator G{x, y), in NLO 
of the l/Af expansion. A thin (green) line represents the classi- 
cal free propagator D{x,y), cf. Eq. <29> . The (red) crosses repre- 
sent classical fields (j}(x), the thick (blue) lines the full propagator 
G(x,y). The thick wiggly lines represent the resummed nonlocal 
interaction A(2;, y; G) recursively defined in Eq. <65> as illustrated 
in Figure|6| Each end of the solid lines as well as the external fields 
carry, as before, a field index and a space-time argument, which are 
summed/integrated over at the vertices. 



(82) 



For the approximation in Eq. \%\\ to be justified we assume 
that bound states of two and more atoms are unlikely to be 
formed, that three-body collisions can be neglected, and that 
the interaction potential V is constant in time. 

We solve the corresponding system of equations in momen- 
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turn space, 

( - ia2,^k^t - M,k{t,p; F)) Fkj(t,t';p) 

dt" ^^{t,t";V,G)Fkj{t",t'-p) 



dt"j:[f,{t,t";p,G)pkjit",t';p), (83) 
- i<J2,ikdt - Mik{t,p; F)^ pkj{t,t';p) 

dt"Y.P^^{t,t"-p-G)pkj(t",t'-p), (84) 



where 



M,,{t,p- F) = 5,, [|- + ^ / Fu{t, t; k) 
+ 51D / F,j{t,t;k). 



(85) 



The non-local self-energies Yi^-''(t, t';p), to NLO in the l/Af 
expansion, are given in AppendixlEI 

With initial values for Fy (0, 0;p), pij{0,0;p), the above 
coupled system of integro-differential equations with first- 
order time derivative yields the time evolution of the corre- 
lation functions, in particular, for i = r, of the momentum 
distribution 



n{t, p) = i (i, t; p) + F22 {t, t;p)-l 
As initial conditions we chose a Gaussian distribution 



n(0,p) 



(86) 



(87) 



Furthermore we chose the initial pair correlation function (I52> 
to vanish, 

0=(*(p,t)#(p,t))e 

= ^(^Fn{t,t;p)- F22{t,t;p)^ +iFi2{t,t;p), (88) 

for t = 0, in accordance with total atom number conservation 
at non-relativistic energies i96ll . As a consequence, 

Fn(0,0;p) = F22(0,0;p) = n(0,p) + 1/2, (89) 
Fi2(0,0;p)=^^2i(0,0;p) = 0. (90) 

Indeed, the condition ( I88> proved to be conserved in time, 
as required by number conservation, such that n{t,p) = 
Fii{t, t;p) — 1/2. The Bose commutation relations, further- 
more, require that 

pn{t,t;p)^P22{t,t;p) = 0, (91) 
-pi2it,t;p) = p2iit,t;p) = l. (92) 

As a specific example we consider an ultracold gas of ^'^Na 
atoms with a line density of ni = 10^ m~^. The dimension- 
less interaction parameter is chosen to be 7 = 7.5 • 10^'* 
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FIG. 8: (Color online) Time evolution of occupation numbers of 
different momentum modes of a uniform one-dimensional gas of 
interacting ^"^Na atoms, for a dimensionless interaction parameter 
7 = 7.5 ■ 10~* and an initially Gaussian distribution, Eq. <87L with 
width a — 1.3 • 10^ m~^, according to the dynamic equations J83> . 
<84t in next-to-leading order in the 1 /M approximation. Shown are 
the relative populations n[t,pi) /n\L on a logarithmic scale as a 
function of evolution time t, where n{t,p) — Fii{t,t;p) — 1/2 
(cf. Eqs. <86t . <89> and the discussion in the text), and ni = 10^ m~^ 
is the total line density. The momenta pi = 2ni /L in the periodic 
box of length L = 32 pm are labelled by their mode numbers i. The 
lattice spacing is =2 pm such that the figure shows the entire 
momentum spectrum. The dashed lines labelled by "HF" show the 
time evolution of the lowest momentum modes i = 0, 1, 2 resulting 
from the Hartree-Fock approximation to the dynamic equations <83> . 
(84). In this approximation, the momentum distribution is constant 
in time! The full NLO dynamics given by the solid lines shows that 
a quasistationary momentum distribution is reached quickly. 



(Cf. Eq. ( I82» . In an experiment such an interaction param- 
eter could be realized, for a given line density and s-wave 
scattering length, by confining the gas in the transverse di- 
rection by a sufficiently steep potential. Assuming a trans- 
verse harmonic potential, the chosen value for 7 corresponds 
to an oscillator length £±_ — -y/2a/ni7 — 0.88 pm, where 
a — 54.6 flBohr = 2.89 nm is the s-wave scattering length of 
the ^"^Na atoms in the commonly trapped state \F ~ l.nip ~ 
— 1) lfl4ll. The width of the initial momentum distribution is 



taken to be cr = 1.3 • 10^ m"\ 

Figure |8l shows the time evolution of the occupation num- 
bers n{t,pi), normalized to the total number of particles in 
the box, niL, according to the dynamic equations (I83t . (I84> 
in next-to-leading order in the 1/J\f approximation. The rel- 
ative populations n{t,pi)/niL are shown on a logarithmic 
scale. The momenta pi = 27ri /L in the periodic box of length 
L — 32 pm are labelled by their mode numbers i. The lattice 
spacing is Os = 2 pm, so the figure shows the entire momen- 
tum spectrum. 

We point out that the dynamic equations ( I83> . (I84> in the 
Hartree-Fock (HF) approximation, where the non-local self- 
energies vanish identically, do not give any change of 
the initial momentum distribution in time. This is due to the 
absence of direct scattering in the HF-approximation. The 
HF- as well as the LO XjM approximation suffer from an 
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FIG. 9: (Color online) Envelopes of the two-point statisti- 
cal and spectral functions Fij{t,0;pi) and Pab{t,o;pi), a,b £ 
{1,2}, respectively, normalized to the total number of particles 
in the box niL, as a function of time t for different momen- 
tum modes pi = 27ri/L. These envelopes have been calcu- 
lated as F{t,0;p,) = [^£l^^F-,a{t,0■,p^fY^'' and p{t,o-p,) = 

Ea=ii?i"(^' ^'iPi)^V^^' respectively. All parameters are chosen as 
in Fig.|8| The statistical functions F are depicted by thick lines, with 
pattern (and color) labelling a particular momentum mode. The spec- 
tral functions p are given by the correspondingly drawn thin lines. 
One observes that at time t — the normalized spectral functions all 
go to the fixed value l/(2niL), in accordance with the equal time 
commutator. At times t > 0.01 s, the respective statistical and spec- 
tral functions for a particular momentum mode show an exponential 
behaviour with the same decay constant 7damp(pi)- 



infinite number of spurious conserved quantities, which are 
not present in the fully interacting theory. These spurious 
constants of motion are associated to an infinite life-time of 
quasi-particle momentum modes, which prevent relaxation to 
a thermal distribution 1 33, 84]. In Figure|8]this is indicated by 
the dashed lines labelled "HF" which give the occupation of 
the lowest three momentum modes i = 0,1,2. 

In contrast to this, the full NLO dynamics given by the solid 
lines shows that a quasistationary momentum distribution is 
reached after a short evolution time |97]. The characteristic 
time Tdamp of the damping of initial oscillations corresponds 
to the time over which the detailed information about the ini- 
tial state is lost. We have checked explicitly that for different 
initial momentum distributions, which are peaked around a fi- 
nite value of \p\, but have a smaller width, such that their total 
particle number and energy content are the same, the momen- 
tum distribution settles to the same quasistationary distribu- 
tion within Tdamp- 

This loss of information shows up in the decrease of the 
two-point functions at large time differences t — t', where 
the time dependence of its envelope becomes an exponential 
of < — Figure |9] shows, for different momentum modes 
Pi = 27Ti/L, the envelopes of the two-point statistical and 
spectral functions Fij{t, Q;pi) and pab{t, 0;pi), a,b E {1, 2}, 
respectively, normalized to the total number riiL, as a func- 
tion of time t. One observes that at times greater than about 
0.01 s, the respective statistical and spectral functions for a 
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FIG. 10: Decay time Tdamp ~ Tdamp derived from the slope 
— 7damp/2 of the envelopes of the two-point statistical and spectral 
functions on the logarithmic scale shown in Fig.|9| for the different 
momentum modes, at time t = 0.02 s. Tdamp is shown as a function 
of momentum, p — 2 sm{asPi/2) / as is the momentum correspond- 
ing to the discretized Laplacian for a lattice constant as- Within 
the resolution of the figure, the results for the statistical and spec- 
tral functions are identical. Tdamp is the characteristic decay time 
over which specific information about the initial state is lost. This 
corresponds to the time when, for a particular momentum mode, the 
strong initial oscillations seen in Fig.[8|are damped out, and the sys- 
tems enters the quasistationary period. The time scale Tdamp is, in the 
weak-coupling regime considered here, orders of magnitude smaller 
than the time scale Teq for the thermalization of the system. 



particular momentum mode show an exponential behaviour 

^exp{-7damp(p)(i-i')/2}, (93) 

with the same decay constant 7dainp(Pi)- Figure^] we 
show the corresponding time scale Tdamp (p) = 7damp(p) 
the different momentum modes pi, at time t = 0.02 s. Within 
the resolution of the figure, no difference between the results 
for the statistical and spectral functions appear. Comparing 
with Fig. [slit becomes clear that Tdamp (p) corresponds to the 
time when, for a particular momentum mode, the strong initial 
oscillations of the occupation number are damped out, and the 
systems enters the quasistationary period. However, the time 
scale Tdamp is, in the weak-coupling regime considered here, 
orders of magnitude smaller than the time scale t^q for the 
thermalization of the system, cf. also Refs. 133113^ ISjI . 

It is important to realize that for times smaller than 
Tdamp (p) the statistical and spectral functions are not char- 
acterized by the same damping. This reflects the fact that 
the far-from-equilibrium dynamics is not characterized by a 
fluctuation-dissipation relation, which would relate the damp- 
ing of F and p (see Eq. i26i). On the other hand, our results 
indicate the approximate validity of a fluctuation-dissipation 
relation for times larger than Tdamp- 

The occupation numbers of the different momentum modes 
are shown in Figure^J for the same parameters as in Fig.|8l 
laid out vs. the momentum p = 2sin{asPi/2)/as correspond- 
ing to the discretized Laplacian for a lattice constant a^. The 
relative populations n{t,pi)/niL are again shown on a loga- 
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FIG. 11: (Color online) Occupation numbers of the modes 
vs momentum, for consecutive time steps of the evolution of 
the one-dimensional gas. Shown are the normalized populations 
n{t,pi) /n\L on a logarithmic scale as a function of momentum p, 
for a box length L = 32 /im and a lattice spacing as = 1.73 pm. 
The (red - color online) circles show the far-from-equilibrium ini- 
tial Gaussian momentum distribution, Eq. j87t . the (black) squares 
the momentum distribution at t = ta = 0.03 s, and the (brown) tri- 
angles the distribution at times = v x 0.01s, i/ — 1,2 (their 
linear interpolation being added as a guide to the eye). The latter, 
quasistationary distribution developed after a short time is close to 
the thermal equilibrium distribution, which is expected, however, to 
be reached only much later. The (blue) diamonds represent an esti- 
mate for the true equilibrium distribution corresponding to the initial 
density and energy. 



rithmic scale such that a Gaussian corresponds to an inverted 
parabola. 

The (red) circles show the far-from-equilibrium initial 
Gaussian momentum distribution, Eq. iEll . the (black) 
squares the momentum distribution at < = ~ 0.03 s (their 
linear interpolation being added as a guide to the eye). The 
(brown) triangles interpolated by dashed lines which corre- 
spond to intermediate times ti = 0.01s and t2 = 0.02 s 
demonstrate that the distribution quickly converges. No sig- 
nificant change of the momentum distributions have been 
found for a larger momentum cutoff l/a^. From the initial 
momentum distribution we have calculated the initial total en- 
ergy 



Pi ki 

^'°^Fa6(0,0;fc,)}^ah(0,0;p,) 



+ 



2 

^ ^"-(O^P) + gmni{niL + 1). 



(94) 



Given this energy and total density ni, as well as the fact that 
the gas is weakly interacting (riia ~ 0.03 ^ 1) we derived 
an estimate for the equilibrium momentum distribution in the 
HFB approximation, ncq(p) — (Wp + Wp){exp[/3((jj(p) — /x)] — 
1}^^ + Wp, where /i is the chemical potential, [3 the inverse 
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FIG. 12: (Color online) The temperature variable Q{tn,p) derived, 
using Eq. <95t . from the non-equilibrium momentum distribution at 
times tu = V X 0.01 s, u = 1, 2, 3 shown in Fig. II II Since a 
quadratic, i.e., free dispersion relation, ijj{p) — p^ /2m, has been as- 
sumed, the data indicates the curvature of inverted parabolae locally 
fitted to the momentum distributions at different times t„ and mo- 
menta Pi. For t — > oo, the values of Q are expected to approach the 
temperature T = 37 nK of the equilibrium distribution indicated by 
the (blue) solid line. 



temperature, Up = {1 - \{v^ 1^2. — u(j>))l^n\ + 1]^} 
and Vp = {[(pV2 - + 1]-^ _ l}-i/2 the HFB 

coefficients, and is defined by the Bogoliubov disper- 
sion relation = {p/m)^J^n\ +p'^ /A. This estimate 
was obtained using number and energy conservation, i.e., 
niL = Y.pi "cqfe) and = J2p, t^(p)"-oq(Pi), respec- 
tively. The resulting equilibrium distribution is shown in Fig. 
II llbv the (blue) diamonds. Hence, as can be seen in the Fig- 
ure, the quasistationary momentum distribution reached at the 
end of the evolution shown in Fig.|8]is already close to the 
expected equilibrium distribution. 

In experiment, a usual procedure for determining the tem- 
perature of an ultracold gas is to fit a Gaussian to the observed 
momentum distribution of the atoms. It is therefore inter- 
esting to compare the temperature one would obtain for the 
non-equilibrium distribution during the quasistationary period 
with the temperature of the anticipated equilibrium distribu- 
tion indicated by the diamonds in Fig. ^2 For sufficiently 
large momenta p the dispersion relation of the weakly inter- 
acting atoms will be close to that of free particles: uj{p) ~ 
p^ /2m. Hence, assuming a Bose-Einstein distribution we de- 
fine an effective temperature variable as 



Q{t,p) = { 



mkB d 
p dp 



In 



1 



7l{t,p) 



(95) 



This is shown, as a function of p, for different evolution times 
= V y. 0.01s, V — 0, ...,3, in Figure [T2l and compared 
to the temperature T of the anticipated equilibrium distribu- 
tion (diamonds in Fig.llU. For large p, where the dispersion 
is approximately quadratic, the value for 6 derived from the 
equilibrium distribution is far from the equilibrium tempera- 
ture T = 37 nK. However, at small momenta, where the oc- 
cupation numbers are large and can more easily measured in 
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experiment, the order of magnitude of the derived temperature 
variable 8(p) coincides with that of the anticipated tempera- 
ture. In contrast to this, the values Q derived from the non- 
equilibrium distribution deviate significantly from the equilib- 
rium temperature over the entire range of momenta and times. 
Hence, our results show, that during the non-equilibrium qua- 
sistationary period which sets in very quickly, the amount of 
information about the final temperature which can be deduced 
from the momentum distribution of the gas is limited. 



VI. CONCLUSIONS 

We have presented a dynamical many-body theory of an ul- 
tracold Bose gas which systematically extends beyond mean- 
field and perturbative quantum-field theoretical procedures. 
The approach is valid for large interaction strengths and/or 
long evolution times, two conditions which have gained a 
great importance in experiments with ultracold atomic gases. 
Examples include the Feshbach-enhanced collisions, dilute 
gases in one- and two-dimensional trapping geometries, the 
Superfluid-to-Mott-insulator transition in optical lattices, and 
the transition from a Bose-Einstein- to a Bardeen-Cooper- 
Schrieffer-type superfluid in a degenerate Fermi gas. The 
methods described in this article are based on an expansion 
of the two-particle (2P1) effective action in powers of the in- 
verse I /J\f of the number of field components. 

To our knowledge, this work represents the first full treat- 
ment of the nonrelativistic dynamics of ultracold atomic gases 
to next-to-leading order in the 2P1 l/7V^-expansion. We have 
reviewed the general framework for obtaining, from the 2P1 
effective action, (exact) dynamic equations for the lowest 
order many-body correlation functions, the mean field and 
the two-point functions which define the non-condensate and 
anomalous density matrices. To be able to solve these equa- 
tions we discussed different approximation schemes, which, 
in the present approach, are implemented at the level of the ef- 
fective action. This ensures the dynamic equations to preserve 
crucial symmetries such as leading to energy and number 
conservation. We first described how the well-known time- 
dependent Hartree-Fock-Bogoliubov (HFB) equations which, 
in the recent past, have been used extensively in the context 
of non-equilibrium dynamics of ultacold atomic gases, result 
as the lowest-order approximation. In the subsequent sections 
we have described the setting up of the non-perturbative 2P1 
dynamic equations to next-to-leading order in the 1/JV ex- 
pansion and their application to the dynamics of an ultracold 
quantum-degenerate gas of sodium atoms in one spatial di- 
mension. As a specific example we have studied the time evo- 
lution of a uniform gas with an initially far-from-equilibrium 
momentum distribution. We have shown that the many-body 
system quickly loses the precise information about the de- 
tails of its initial state and enters a quasi-stationary evolu- 
tion period. During this period, the momentum distribution 
is still far from equilibrium, and the drift towards the antici- 
pated Bose-Einstein equilibrium distribution is found to be ex- 
tremely slow. Fitting a Bose-Einstein distribution during the 
quasi-stationary period yields effective temperatures far away 



from the expected final temperature. 

Summarizing, we point out, that, in order to describe non- 
equilibrium dynamics of strongly interacting and/or dense 
atomic gases, an approach beyond the conventional Hartree- 
Fock(-Bogoliubov) approximation is required. In fact, in the 
case where the system is described solely by the two-point 
functions, i.e., whenever the mean field vanishes, the HF dy- 
namic equations do not describe any dynamical change in the 
mode occupation numbers such that thermal equilibration is 
not accounted for This is due to the fact, that the HF approx- 
imation neglects multiple scattering and therefore conserves 
an infinite number of spurious conserved quantities. 

We emphasize that the 2P1 1/J\f expansion to next- to- 
leading order represents an approximation scheme which is 
capable of describing the dynamics of strongly interacting 
many-body systems far from thermal equilibrium. The study 
of its implications in the context of ultracold atomic gases is 
an exciting and demanding task for near-future research. 
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APPENDIX A: DERIVATION OF T [<?!), G] TO ONE-LOOP 
ORDER 

In Eq. jl9> the 2PI effective action is explicitly given up to 
one-loop order while all higher order contributions are de- 
noted by r2[0, G]. The one-loop expression can be derived 
as follows: We first note that the 2PI effective action can be 
viewed as the IPl effective action 

r^[(/.] - W[J,K] - j ct>,{x).h{x) (Al) 

J X 

for a theory governed by the modified classical action 

= 5[<i>] + i / ^,{x)K,,{x,y)<^,{y), (A2) 

^ J xy 

which contains an additional "potential" term quadratic in the 
fields. For K = Q, [(f)] is equivalent to the conventional IPI 
effective action. The IPl effective action r^[0] then assumes, 
to one-loop order, the well-known form 

where iGQlj{x, y; (f)) — 5"^ S[(l)]/ 5(l)i{x)5(l)j{y) is the classical 
inverse propagator 
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Using the definition iAl\ one finds that 

Sr^[(j)] 6W[J,K] 



(A4) 



and, together with Eq. ilbi . the Legendre transform of [(f)] 
with respect to K, i.e., the 2PI effective action can be written 
as 



<5r^ [(/)] 



■Rji{y,x) 



r'^M-^ / Mx)K,,{x,y)^,{y) 

^ J xy 



- -TtKG. 
2 



(A5) 



Now we insert Eq. ( IA3> into Eq. ( IA5> to obtain the 2PI effec- 
tive action in one-loop approximation: 

r"°°P[c6. G] = S\<^\ + iTr(i In Gg ^[0] - - /^g) . 

(A6) 



exact inverse propagator G ^ and the classical inverse propa- 
gator of the modified action , 



iGijix,y) 



(5^r^ [0] 
S(f>i{x)S(f>j{y) 

Goh ix,y)- iKij {x, y) - Hf^ (x, y) 



(A7) 

with Tifj {x,y) being the proper self energy to which only one- 
particle irreducible diagrams contribute, one obtains, from 
Eq. ( IA6> . the one-loop terms of r[(j>, G] in Eq. ( I19> . 



APPENDIX B: THE FUNCTIONS Mij, T,^°\ Eg, AND Ef^ IN 
NLO OF THE 1 /Af EXPANSION 



In Section IIVBI we have derived the dynamic equations for 
the real fields 4>i{x), Fij{x,y), and pij{x,y) in NLO of the 
1/JV expansion of the 2P1 effective action, cf. Eqs. ( I78t . 
il9i . and ( I80l l, respectively. In these equations, the real 
functions M,j (x; </>, F), E^f {x, y; 0, F), Eg (x, y; G), and 



Using the one-loop relation G ^ = Gq ^ — between the E^ (a;, 0, G) are given in terms of 0, F, and p as follows: 



(a:,z/;(/), i^) = (5ij(5c(a; - y) iJiB(a;) + ^ J V{x - z)(^(f>k{z)(j)k{z) + Fkk{z, z)^ + V{x - y)(j)i{x)(j)j{y) + Fij{x,y)^ 



(Bl) 



Sg(a;,j/;(^,G) = -AF{x,y) (j)i{x)(j)j{y) + F,j{x,y) + -Ap{x,y)pij{x,y) - Pp{x,y) Fij{x,y) + -Pp{x,y)) pij{x,y), 



^tj{x,y;(l),G) = -Ap{x,y) (l)^{x)(t)j{y) + Fij{x,y) ~ Ai?(x, y) (x, y) - Pp{x,y) F^j{x,y) - Pp {x , y)) p^j {x , y) . 



(B2) 
(B3) 



In here, the functions Af,p{x, y) are given as 



AF,p(a;, y) ^ Jf J If,p{x, z) V{z - y), 
If{x, y)^J^ V{x - z) (f{z, yf - yf) 



dz I, 



,{x,u)V{u-z) (i^(z,y)2- lp(z,y)2) 



dz If{x, u) V{u - z) Fij {z, y)pij (z, y) 



-fp(a;,y) = 2| f V{x - z)Fij{z,y)p,j{z,y) - f f dz Ip{x,u)V{u ~ z) F,j{z,y)pij{z,y)^ 



(B4) 

(B5) 
(B6) 
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where F"^ = FijFij, etc., and the functions Pp^p{x, y) as 

Pp{x,v) = ^ f { f V{x-v)HF{v,w)V{w-y) 



' V J w 

rva 



dw 



Xo 



dw 



V{x - v) Hf{v,w) Ip{w, y) + If{x, w) Hp{w, v) V{v - y) 
V{x ~ v) Hp{v,w) lF{w,y) + Ip{x,w) Hf{w,v)V{v - y) | 



10 

dv / dw Ip{x,v) Hf{v,w) Iplwjy) + / dv dw Ip{x,v) Hp{v,w) iFiWjy) 
Jo Jo Jo Jo 

rVo rvo 

+ dv dw If{x,v) Hp{v,w) Ip(w,y), 

Jo J zo 

4 



(B7) 



V{x — v) Hp{v, w) V{w — y) 



f V >J w 

(•xo 



dw 



V{x - v) Hp{v,w) Ip{w,y) + Ip{x,w) Hp{w,v) V{v - | 



Vo 

Xo nvo 

dv dw Ip(x,v) Hp{v,w) Ip{w,y) 



(B8) 



with 



HF{x,y) = (l),{x)F,,{x,y)4>j{y), (B9) 
Hp{x,y) = (t)i{x) p^j{x,y)(l)j{y), (BIO) 

where it is, as usual, summed over double indices. In 
Eqs. iBll and ( IB8> it has been used that all time integrals over 
the closed path C, of HfIf, IfHf, and IfHpIf vanish. 

Finally, we provide the momentum-space self energies 
'E,fj''^{t,t';p) in 1+1 dimensions which enter the dynamical 
equations ( I83t . ( I84t . From Eqs. (IB2t to ( IB6t . one obtains, for 
a homogeneous system, with (pi = 0, and a coupling given in 
Eq. (18 1> . by Fourier transformation: 



i:[At,t';p;G) = -giB / lF{t,t' -p ~ k)F,,{t,t' -k) 



-Ip{t,t']p- k) p,j{t,t';k) 



(BID 



Ef,(t,t';p;G) --.giD / Ip{t,t' ■,p ~ k)F,,{t,t' ■ k) 



+ lF{t,t';p- k)pij{t,t';k) , (B12) 



lF{t,t';p)=giB / \F,j{t,t';p-k)F,j{t,t';k) 



dt" Ip{t,t":p- k) 



X (F,,{t\t'-k-k')F,,{t",t'-k') 
-^p,,{t",t';k-k')p,j{t",t';k')) 



2 / dt" lF{t,t"-p-k) 



X K 



,{t,t'-k-k')p,j{t,t'-k')\Y (B14) 



Here, J^, = (27r) ^ J dk denotes the one-dimensional mo- 
mentum integral. 



Ip{t,t';p)^2giD / \F,j{t,t';p-k)p,j{t,t';k) 



k' Jf 



dt" Ip{t,t";p- k) 

X F,,{t,t'-k~k')p,,{t,t'-k')Y (B13) 
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[93] In Refs. |6^ j^ll, (sj, the following notation was used: 
*(x;t) = </>(x,i), r(x,y;t) = n(x, y; t), and $(x, y; t) = 
7n(x,y;f). 

[94] For example, N = 2 gives 5 independent invariants, in accor- 
dance with the two field components of (j) and the three inde- 
pendent components of the Hermitian matrix G. 

[95] We divide the self-energy into a 'local' and a 'non-local' part, 
E = E*"' + S. For a local coupfing V{x - y) = VS{x - y), 
the contribution E'°' , which we include in M, would be a local 
function, E™ {x, y) = E'°' {x)5{x - y) Q. 

[96] At non-relativistic energies, the total number of atoms in a cer- 
tain finite system A is fixed, say it is TV. Then, for any subsys- 
tem B of A, the density matrix, due to number conservation, 
takes the form p = Eiv^.iv^ \N~Nb){N-N'b\®\Nb){N'b\. 
Hence, like the mean field (^'), the expectation value (^'^') = 
Tr Tr necessarily vanishes, where Tr a\b de- 

notes the trace over all basis states of the Hilbert space cor- 
responding to the "environment" of B in A. 

[97] Other than in a relativistic theory, where 1-^3 on-shell re- 
actions with pair production lead to thermalization js^, at 
non-relativistic energies strict number conservation only admits 
2 — » 2- vertices. Moreover, due to the restricted phase space in 
1-1-1 dimensions, momentum conservation leaves all momenta 
invariant. Hence, the momentum occupations are not damped 
in 2 — > 2 on-shell reactions, and thermalization requires a finite 
width of the propagator. 



